#include <iostream>
#include <fstream>
#include <vector>
#include <cmath>
struct DataPoint {
    double x, y;
};

// 读取数据文件
void read_data(const std::string& filename, std::vector<DataPoint>& data) {
    std::ifstream in_file(filename);
    if (!in_file) {
        std::cerr << "Error opening file!" << std::endl;
        return;
    }

    double x, y;
    while (in_file >> x >> y) {
        data.push_back({ x, y });
    }

    in_file.close();
}

// 计算预测值
double predict(double x, const std::vector<double>& params, int t) {
    if (t == 1) {  // 一次函数
        return params[0] * x + params[1];
    }
    else if (t == 2) {  // 二次函数
        return params[0] * x * x + params[1] * x + params[2];
    }
    else if (t == 3) {  // 三次函数
        return params[0] * x * x * x + params[1] * x * x + params[2] * x + params[3];
    }
    return 0.0;
}

// 计算梯度下降的损失函数和梯度
void gradient_descent(const std::vector<DataPoint>& data, std::vector<double>& params, int t, double learning_rate = 0.001, int iterations = 1000) 
{
    int m = data.size();
    const double epsilon = 1e-8;  // 避免除以零
    const double gradient_clip = 100.0;  // 防止梯度爆炸
    for (int iter = 0; iter < iterations; ++iter) 
    {
        // 初始化梯度
        std::vector<double> gradients(params.size(), 0.0);

        // 计算梯度
        for (const auto& point : data) 
        {
            double error = predict(point.x, params, t) - point.y;
            if (t == 1) 
            {  // 一次函数
                gradients[0] += error * point.x;  // 对a求偏导数
                gradients[1] += error;            // 对b求偏导数
            }
            else if (t == 2)
            {  // 二次函数
                gradients[0] += error * point.x * point.x;  // 对a求偏导数
                gradients[1] += error * point.x;            // 对b求偏导数
                gradients[2] += error;                      // 对c求偏导数
            }
            else if (t == 3) 
            {  // 三次函数
                gradients[0] += error * point.x * point.x * point.x;  // 对a求偏导数
                gradients[1] += error * point.x * point.x;            // 对b求偏导数
                gradients[2] += error * point.x;                      // 对c求偏导数
                gradients[3] += error;                                // 对d求偏导数
            }
        }

        // 更新参数
        for (int i = 0; i < params.size(); ++i) 
        {
            params[i] -= learning_rate * gradients[i] / m;
        
        // 进行梯度裁剪
            if (std::abs(gradients[i]) > gradient_clip) {
                gradients[i] = std::copysign(gradient_clip, gradients[i]);
            }
        }

    // 检查是否出现 NaN
        for (double param : params) 
        {
            if (std::isnan(param)) 
            {
                std::cerr << "NaN encountered in parameters during iteration " << iter << std::endl;
                return;
            }
        }
    }
}

// 计算残差平方和（拟合度）
double compute_residual_sum(const std::vector<DataPoint>& data, const std::vector<double>& params, int t) {
    double sum = 0.0;
    for (const auto& point : data) {
        double error = predict(point.x, params, t) - point.y;
        sum += error * error;
    }
    return sum;
}